clear all
use  "[Folder]\Panel_Firm_Lab_DirVarOcc_temp.dta" 

	global MLab = 11
	global FLab = 4
	distinct Ind if u==0
	global Inds = r(ndistinct) - 1
	global Sizes = 4
	global Cities = 42
	global Years = 6
	global Res = 14
	global NbrOfParam =  $MLab + $FLab + $Cities + $Sizes + $Years + $Inds + $Res
	disp $NbrOfParam


//Function evaluation program
     program nlceslabrho
                version 14.2
                syntax varlist [if], at(name)
                local log_VA: word 1 of `varlist'
				
	// Retrieve parameters out of at matrix				
				tempvar LabAggH LabAggL  temp_VA LabAggH_E LabAggL_E LabAggH_C LabAggL_C
				gen double `LabAggH' = 0 `if'
				gen double `LabAggL' = 0 `if'

				gen double `LabAggH_E' = 0 `if'
				gen double `LabAggL_E' = 0 `if'
				gen double `LabAggH_C' = 0 `if'
				gen double `LabAggL_C' = 0 `if'

                gen double `temp_VA' = 0 `if'
				
				replace `LabAggL_E' = `LabAggL_E' + RefGr  `if'

				//Male lab types
				forvalues i=1(1)3 {
				tempname delta 
				scalar `delta' = `at'[1,`i']
				replace `LabAggH_E' = `LabAggH_E' + `delta' * (lab_MH_`i')  `if'
				}
				
				forvalues i=1(1)2 {
				tempname delta 
				scalar `delta' = `at'[1,3+`i']
				replace `LabAggL_E' = `LabAggL_E' + `delta' * (lab_ML_`i')  `if'
				}
				
				forvalues i=1(1)3 {
				tempname delta 
				scalar `delta' = `at'[1,3+2+`i']
				replace `LabAggH_C' = `LabAggH_C' + `delta' * lab_MH_C`i'  `if'
				}

				forvalues i=1(1)3 {
				tempname delta 
				scalar `delta' = `at'[1,3+2+3+`i']
				replace `LabAggL_C' = `LabAggL_C' + (`delta' * lab_ML_C`i')   `if'
				}
				//Female lab types
				
				tempname  delta 
				scalar `delta' = `at'[1,$MLab +1]
				replace `LabAggH_E' = `LabAggH_E' + `delta' * (lab_FNat_H) `if'
				tempname  delta 		
				scalar `delta' = `at'[1,$MLab +2]
				replace `LabAggH_E' = `LabAggH_E' + `delta' * (lab_FFor_H)  `if'

				tempname  delta 
				scalar `delta' = `at'[1,$MLab +3]
				replace `LabAggL_E' = `LabAggL_E' + `delta' * (lab_FNat_L ) `if'
				tempname  delta 
				scalar `delta' = `at'[1,$MLab +4]
				replace `LabAggL_E' = `LabAggL_E' + `delta' * (lab_FFor_L )  `if'
				

				
				
				//Adj index for number of industry dummies
				forvalues i=1(1)$Inds {
				tempname Industry 
				scalar `Industry' = `at'[1,$MLab + $FLab + `i']
				replace `temp_VA' = `temp_VA' + `Industry'*Ind`i' 	`if'			
				}

				//Adj index for number of firm size dummies
				forvalues i=1(1)$Sizes {
				tempname Size 
				scalar `Size' = `at'[1, $MLab + $FLab  + $Inds +`i']
				replace `temp_VA' = `temp_VA' + `Size'*FirmSize_`i' `if'				
				}
				//Adj index for number of City dummies
				forvalues i=1(1)$Cities {
				tempname City 
				scalar `City' = `at'[1, $MLab + $FLab  + $Inds + $Sizes +`i']
				replace `temp_VA' = `temp_VA' + `City'*k_`i' 	`if'			
				}
				//Adj index for number of year dummies
				forvalues i=1(1)$Years {
				tempname Year 
				scalar `Year' = `at'[1, $MLab + $FLab  + $Inds + $Sizes +$Cities +`i']
				replace `temp_VA' = `temp_VA' + `Year'*y_`i' `if'
				}				
				
		
				
				
				tempname rho_H
				scalar `rho_H' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +13]
				replace `LabAggH' = `LabAggH_E'^(`rho_H') + `LabAggH_C'^(`rho_H')   `if'
				
				tempname rho_L
				scalar `rho_L' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +14]
				replace `LabAggL' = `LabAggL_E'^(`rho_L') + `LabAggL_C'^(`rho_L')   `if'
		

				
				tempname  ic thetaL Cap Cap_Sq Cap_Qu TI TI_Sq TI_Qu TICap TICap_Sq TI_SqCap TI_SqCap_Sq   rho_H rho_L 
				scalar `thetaL' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +1]
				scalar `Cap' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +2]
				scalar `Cap_Sq' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +3]
				scalar `Cap_Qu' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +4]
				
				scalar `TI' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +5]
				scalar `TI_Sq' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +6]
				scalar `TI_Qu' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +7]
				
				scalar `TICap' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +8]
				scalar `TICap_Sq' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +9]
				scalar `TI_SqCap' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +10]
				scalar `TI_SqCap_Sq' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +11]
				scalar `ic' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +12]

				scalar `rho_H' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +13]		
				scalar `rho_L' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +14]	
				
				replace `temp_VA' = `temp_VA' + `ic'  + `thetaL' *ln( `LabAggH'^(1/`rho_H') + `LabAggL'^(1/`rho_L') )  + `Cap'*logCap + `Cap_Sq'*logCap_Sq + `Cap_Qu'*logCap_Qu + `TI'*logGrosInvTot + `TI_Sq'*logGrosInvTot_Sq  + `TI_Qu'*logGrosInvTot_Qu + `TICap'*logGrosInvTot*logCap + `TICap_Sq'*logGrosInvTot*logCap_Sq  + `TI_SqCap'*logGrosInvTot_Sq*logCap + `TI_SqCap_Sq'*logGrosInvTot_Sq*logCap_Sq  `if'
				replace `log_VA' = `temp_VA' `if'
            end
			///////////////////PROGRAM END///////////////////////////
			




			
//run regression, with number of param matching M			
mat M = J(1,$NbrOfParam, 0.7)
mat M[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +12]=5



nl ceslabrho @ logFtg_Foradlingsvarde if u==0, nparameters($NbrOfParam) initial(M) iter(100) variables(lab_MH_C1-lab_MH_C3 lab_ML_C1-lab_ML_C3 lab_ML_1-lab_ML_2 lab_MH_1-lab_MH_3 lab_FNat_H lab_FFor_H lab_FNat_L lab_FFor_L  logCap logCap_Sq logCap_Qu logGrosInvTot logGrosInvTot_Sq  logGrosInvTot_Qu logGrosInvTot y_1-y_6 k_1-k_41 FirmSize_1-FirmSize_4  Ind* RefGr)


margins, dydx(RefGr lab_MH_C* lab_ML_C*) atmeans

////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////
//Parameter values in table 8 are marginal effects for respective labor type normalized by the marginal effect of the reference group (except for rho_X which is the parameter value from the original estimation)
////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////
